library(janitor)
Attaching package: ‘janitor’
The following objects are masked from ‘package:stats’:
chisq.test, fisher.test
library(readODS)
library(httr)
library(here)
here() starts at /Users/tomdavie/Documents/GitHub/ev_climate_change_project
# loading in shape file
uk_shape_file <- st_read(here("raw_data/LA_shape/Local_Authority_Districts__April_2019__UK_BFE_v2.shp")) %>%
clean_names() %>%
st_simplify(dTolerance = 1000) %>%
st_transform("+proj=longlat +datum=WGS84") %>%
dplyr::select(lad19cd, long, lat, geometry)
Reading layer `Local_Authority_Districts__April_2019__UK_BFE_v2' from data source
`/Users/tomdavie/Documents/GitHub/ev_climate_change_project/raw_data/LA_shape/Local_Authority_Districts__April_2019__UK_BFE_v2.shp'
using driver `ESRI Shapefile'
Simple feature collection with 382 features and 10 fields
Geometry type: MULTIPOLYGON
Dimension: XY
Bounding box: xmin: -116.1928 ymin: 5333.603 xmax: 655989 ymax: 1220302
Projected CRS: OSGB 1936 / British National Grid
# Load in clean Electric Vehicles by Local Authority data
uk_ev <- read_csv(here("clean_data/ev_by_la_clean.csv"))
Rows: 475 Columns: 40
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (40): ons_la_code_apr_2019, region_local_authority_apr_2019_3, x2021_q1, x2020_q4, x2020_q3, x2020_q2, x2020_q1, x2019_q4, x2019_q...
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
# Joining uk_ev + shape file
uk_ev_map <- uk_ev %>%
left_join(uk_shape_file, by = c("ons_la_code_apr_2019" = "lad19cd")) %>%
drop_na() %>%
mutate(across(c(x2021_q1:x2011_q4), as.numeric)) %>%
st_as_sf()
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
Warning: Problem with `mutate()` input `..1`.
ℹ `..1 = across(c(x2021_q1:x2011_q4), as.numeric)`.
ℹ NAs introduced by coercion
# Set colours and bins
pal <- colorBin("Greens", domain = uk_ev_map$x2021_q1, bins = c(0, 500, 1000, 2500, 5000, 10000, 15000))
# Set labels
uk_ev_map_labels <- sprintf(
"<strong>%s</strong><br/>%g Electric Vehicles",
uk_ev_map$region_local_authority_apr_2019_3, uk_ev_map$x2021_q1) %>%
lapply(htmltools::HTML)
# Geospatial of EV Vehicles in the UK 2021 Q1
uk_ev_map %>%
leaflet() %>%
setView(lng = -4.2026, lat = 55.8, zoom = 4.7, options = list()) %>%
addProviderTiles(providers$CartoDB.Positron) %>%
addPolygons(fillColor = ~pal(x2021_q1),
weight = 0.1,
opacity = 0.9,
color = "black",
fillOpacity = 0.8,
highlightOptions = highlightOptions(color = "green", weight = 2,
bringToFront = TRUE),
label = uk_ev_map_labels,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = pal, values = ~x2021_q1, opacity = 0.7, title = NULL,
position = "bottomright")
pal <- colorBin("Greens", domain = uk_ev_map$x2021_q1, bins = c(0, 500, 1000, 2500, 5000, 10000, 15000))
uk_ev_map_labels <- sprintf(
"<strong>%s</strong><br/>%g Electric Vehicles",
uk_ev_map$region_local_authority_apr_2019_3, uk_ev_map$x2021_q1) %>%
lapply(htmltools::HTML)
# Geospatial of EV Vehicles in the UK 2021 Q1
uk_ev_map %>%
leaflet() %>%
setView(lng = -4.2026, lat = 55.8, zoom = 4.7, options = list()) %>%
addProviderTiles(providers$CartoDB.Positron) %>%
addPolygons(fillColor = ~pal(x2021_q1),
weight = 0.1,
opacity = 0.9,
color = "black",
fillOpacity = 0.8,
highlightOptions = highlightOptions(color = "green", weight = 2,
bringToFront = TRUE),
label = uk_ev_map_labels,
labelOptions = labelOptions(
style = list("font-weight" = "normal", padding = "3px 8px"),
textsize = "15px",
direction = "auto")) %>%
addLegend(pal = pal, values = ~x2021_q1, opacity = 0.7, title = NULL,
position = "bottomright")
# Wrangling to create an EV count over time plot
uk_ev_longer <- uk_ev %>%
# Pivot longer to get year and count columns
pivot_longer(cols = c(x2021_q1:x2011_q4), names_to = c("year"), values_to = "no_of_ev") %>%
# Filter so we only have UK as a whole data AND we only want final numbers of the year so Q4
filter(region_local_authority_apr_2019_3 == "United Kingdom" & str_detect(year, "q4")) %>%
# Simplify to just show year
mutate(year = case_when(str_detect(year, "2021") ~ "2021",
str_detect(year, "2020") ~ "2020",
str_detect(year, "2019") ~ "2019",
str_detect(year, "2018") ~ "2018",
str_detect(year, "2017") ~ "2017",
str_detect(year, "2016") ~ "2016",
str_detect(year, "2015") ~ "2015",
str_detect(year, "2014") ~ "2014",
str_detect(year, "2013") ~ "2013",
str_detect(year, "2012") ~ "2012",
str_detect(year, "2011") ~ "2011"),
year = as.numeric(year),
no_of_ev = as.numeric(no_of_ev))
uk_ev_longer %>%
ggplot() +
aes(x = year, y = no_of_ev) +
geom_col(fill = "#08a44c") +
scale_x_continuous(breaks = c(2011:2020)) +
scale_y_continuous(breaks = seq(0, 220000, by = 20000), limits = c(0, 220000)) +
labs(title = "\nNumber of Electric Vehicles over time in the UK\n",
x = "\nYear\n",
y = "\nNumber of Electric Vehicles\n") +
theme_minimal()

Row binding grid NO2 data
no2_2010 <- read_csv(here("raw_data/no2_by_grid_2010.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2011 <- read_csv(here("raw_data/no2_by_grid_2011.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2012 <- read_csv(here("raw_data/no2_by_grid_2012.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2013 <- read_csv(here("raw_data/no2_by_grid_2013.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2014 <- read_csv(here("raw_data/no2_by_grid_2014.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2015 <- read_csv(here("raw_data/no2_by_grid_2015.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2016 <- read_csv(here("raw_data/no2_by_grid_2016.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2017 <- read_csv(here("raw_data/no2_by_grid_2017.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2018 <- read_csv(here("raw_data/no2_by_grid_2018.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
no2_2019 <- read_csv(here("raw_data/no2_by_grid_2019.csv"), skip = 6,
col_names = c("uk_grid_code", "x", "y", "no2"))
Rows: 281802 Columns: 4
── Column specification ──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (1): no2
dbl (3): uk_grid_code, x, y
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
# Create empty data frame
no2_all <- data_frame()
Warning: `data_frame()` was deprecated in tibble 1.1.0.
Please use `tibble()` instead.
This warning is displayed once every 8 hours.
Call `lifecycle::last_warnings()` to see where this warning was generated.
# For each year, bind rows to one dataset
for (i in 2010:2019) {
df_name <- paste0("no2_", i)
df_input <- as.name(df_name)
df <- eval(df_input) %>%
mutate(year = i)
no2_all <- bind_rows(no2_all, df)
}
# Remove missing NO2 grid values
no2_clean <- no2_all %>%
filter(no2 != "MISSING")
library(proj4)
proj4string <- "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs"
no2_clean_row <- no2_clean %>%
rowid_to_column()
# Source data
xy <- no2_clean_row %>%
select(x, y, rowid)
# Transformed data
pj <- project(xy, proj4string, inverse=TRUE)
Warning in project(xy, proj4string, inverse = TRUE) :
more than two dimensions found, using first two
latlon <- data.frame(xy, lat=pj$y, lon=pj$x)
final <- merge(no2_clean_row, latlon, by.x = "rowid", by.y = "rowid") %>%
filter(year == 2019) %>%
select(lat, lon, no2)
final <- final %>%
mutate(no2 = as.numeric(no2))
bins <- c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10)
pal <- colorBin("Spectral", domain = final$no2, bins = bins, na.color = "transparent", reverse = TRUE)
leaflet() %>%
addProviderTiles("CartoDB.Positron", options = providerTileOptions(noWrap = TRUE)) %>%
setView(lng = -4.2026, lat = 55.8, zoom = 4.7, options = list()) %>%
addHeatmap(data = final,
lng = ~lon,
lat = ~lat,
intensity = ~no2,
minOpacity = 0.1,
max = 45,
radius = 1,
blur = 1) %>%
addLegend(pal = pal, values = final$no2,
title="Average NO2 Conc")